#!/usr/bin/python
#+
#
#  REQUIREMENTS:
#	numpy >=1.4
#
#  HISTORY:
#     - 30/05/2012 : GOX (MPE). First frozen version. Tested and working. 
#                          Many parameters are hard-coded/ or limited flexibility:
#			  e.g. # of subaps per SH pupils & per DSP cannot be changed.
#     - 18/09/2012 : GOX (MPE) - SVN repository created on 
#				https://code.google.com/p/bcu
#
#-
import matplotlib.pyplot as plt
import numpy as np
import warnings

class SubapDef(object):
	'''
	AIM:
		This Class provides the following tables for the ARGOS BCU (for each DSP):
		- dsc_SubapPixel
		- dsc_NumSubapsReady
		- sc_RemapSlopeVector
		
		In addition to the tables for the BCU, 
		it does provides a 2D array, representing the SH pupil definition and 
		subapertures numbering.

	FUNCTION:
		shiftSHPupils: 
		          allowed to shift ([x,y]) the individual SH pupil and 
		          modify the 3 end-vectors accordingly.
			  The pupils can now be shifted safely within their half CCD 
				(as long as no overlap between pupils)
			  
			  Re-definition of the numbers of subaperture can only be a switch from one to the other 					half of CCD (switch between 257 or 271 subap, incl. extra subap)
	'''
	def __init__(self):
		'''
		 For simplicity/'proven-to-work', a few parameters are here 'hard-coded'.
		'''
		self._sc_ReplyVectorPtr = 524288

		self.nsubap_dsp0=257
		self.nsubap_dsp1=271 # nb of -real- subap
		self.max_nb_subapPerDsp=271 # maximum number of -real- subaperture/dsp

		self.subapSize=8

		self._SH1=np.array((264,264))
		self._SH2=np.array((264,264))
		self._SH3=np.array((264,264))

		self.vectorsBCU  = self.setDefaultVectorsBCU()
		self.arrayHuman  = self.setDefaultArrayHuman()
		print "dsc_SubapPixel vector set"
		self.numSubapReady =self.numSubapToDo()
		self.remapVectorBCU = np.array(self.setRemapVectorBCU(),'int')
		self.remapVector = self.setRemapVector()		
		print "dsc_NumSubapsReady and sc_RemapSlopeVector set"



	def setDefaultVectorsBCU(self):
		option=2
		if option==1:
			self.nsubap_dsp0=257
			self.nsubap_dsp1=271
			centers=[[66+4,66],[66+4,198],[198-4,132+4]]
			#centers=[[60+7+8,60+3+4],[60+7+8,60+3+4+60+10+60],[190,132+4]]
			radii=[15,15,15]
		elif option==2:
			self.nsubap_dsp0=271
			self.nsubap_dsp1=257
			centers=[[60+7+8,60+3+4],[60+7+8,60+3+4+60+10+60],[190,132+4-8]]
			radii=[15,15,15]
		tmp=self._createSubapDef(centers,radii)
	

		return self._array2subap_pix_ptr(tmp)
		
	def setDefaultArrayHuman(self):
		arrayHuman=self._subap_pix_ptr2array(self.vectorsBCU[0],self.vectorsBCU[1])
		return arrayHuman
	
	def getVectorsBCU(self):
		return self.vectorsBCU
	
	def getArrayHuman(self):
		return self.arrayHuman
	
	
	def shiftSHPupils(self, s1,s2,s3):
		''' 
		s1,s2,s3 are the [x,y] shift of each SH pupil images
		n.b.: in python x dimension is the vertical one....
		'''
		_isShiftOk=0
		if np.mod(s3[1],8) != 0:
			warnings.warn("Pupil cannot be shifted by shift != of +-8pixel")
			print "Subap-aperture should be within one Half-pnCCD"
		elif np.abs(s3[1]/8.) >1:
			warnings.warn("Pupil cannot be shifted by more than 8px x-axis")
			print "Shift too high!"
		elif ((s3[1]/8. == 1) & (self.nsubap_dsp1>self.nsubap_dsp0)):
			warnings.warn("Shift too high: you cannot shift in that direction")
			print "Shift too high!"
		elif ((s3[1]/8. == -1) & (self.nsubap_dsp1<self.nsubap_dsp0)):
			warnings.warn("Shift too high: you cannot shift in that direction")
			print "Shift too high!"
		elif ((s3[1]/8. == 1) & (self.nsubap_dsp1<self.nsubap_dsp0)):
			tmp0 = self.nsubap_dsp0
			tmp1 = self.nsubap_dsp1
			self.nsubap_dsp1 = tmp0
			self.nsubap_dsp0 = tmp1
			_isShiftOk=1

		elif ((s3[1]/8. == -1) & (self.nsubap_dsp1>self.nsubap_dsp0)):
			tmp0 = self.nsubap_dsp0
			tmp1 = self.nsubap_dsp1
			self.nsubap_dsp1 = tmp0
			self.nsubap_dsp0 = tmp1
			_isShiftOk=1
		
		if _isShiftOk==1 :
			array=self.arrayHuman
			[xcoord_sh1,ycoord_sh1]=np.where(self._SH1 == 1)
			[xcoord_sh2,ycoord_sh2]=np.where(self._SH2 == 1)
			[xcoord_sh3,ycoord_sh3]=np.where(self._SH3 == 1)
			
			[xcoord_sh1,ycoord_sh1]=[xcoord_sh1 + s1[0],ycoord_sh1 + s1[1]]
			[xcoord_sh2,ycoord_sh2]=[xcoord_sh2 + s2[0],ycoord_sh2 + s2[1]]
			[xcoord_sh3,ycoord_sh3]=[xcoord_sh3 + s3[0],ycoord_sh3 + s3[1]]
			
			array=array*0
			self._SH1*=0
			self._SH2*=0		
			self._SH3*=0
			self._SH1[xcoord_sh1,ycoord_sh1]=1
			self._SH2[xcoord_sh2,ycoord_sh2]=1
			self._SH3[xcoord_sh3,ycoord_sh3]=1
			
			array=self._SH1+self._SH2+self._SH3
			if np.where(array > 1)[1].size !=0:
				warnings.warn("Wrong pupil definitions: Shack-Hartmann are overlapping")
				print "Wrong pupil definitions: Shack-Hartmann are overlapping"
			
			self.vectorsBCU=self._array2subap_pix_ptr(array)
			self.arrayHuman=array
			self.numSubapReady =self.numSubapToDo()
			self.remapVectorBCU = np.array(self.setRemapVectorBCU(),'int')
			self.remapVector = self.setRemapVector()

	def setRemapVector(self):
		'''
		likely not necessary anymore. For human representation.
		'''
		slopeRemapVector = np.array(self.remapVectorBCU,'int')
		slopeRemapVector=slopeRemapVector - np.min(slopeRemapVector) + 4
		header_size=8
		size0=np.size(self.vectorsBCU[0])/64 
		size1=np.size(self.vectorsBCU[1])/64 
		print size0,size1

		if size0 < size1 :
			slopeRemapVector[slopeRemapVector[:] > size0*2-1]  = \
			    slopeRemapVector[slopeRemapVector[:] > size0*2-1] -\
			    (size1-size0)*2 - \
			    (header_size)
                else :
			slopeRemapVector[slopeRemapVector > size1*2-1]  = \
			    slopeRemapVector[slopeRemapVector > size1*2-1] -\
			    (header_size)
                

		remapXslopes= (slopeRemapVector[0::2])/2
		remapYslopes= (slopeRemapVector[1::2]-1)/2
	        return [remapXslopes, remapYslopes]

	def setRemapVectorBCU(self):
		'''
		creates the 'sc_RemapSlopeVector' for the BCU, i.e. the final slope remapping vector.
		a few things are hard-coded and rely on a frozen 'sc_ReplyVector'
		'''
		assert self._sc_ReplyVectorPtr == 524288
		len_hdr_RV=8
		
		subapList_dsp0=self._subap_Numbering(self.arrayHuman[:,0:132]*1,\
							list([[],[],[]]))
		subapList_dsp1=self._subap_Numbering(np.fliplr(self.arrayHuman[:,132:264])*1,\
							list([[],[],[]]))


		subapList_dsp1[0].reverse()
		subapList_dsp1[1].reverse()
		subapList_dsp1[2].reverse()
	
		remap_list=[]

		subapList_dsp0=np.array(subapList_dsp0)
		subapList_dsp1=np.array(subapList_dsp1)

		toto_dsp0=subapList_dsp0[0]*264 + subapList_dsp0[1]
		toto_dsp1=subapList_dsp1[0]*264 + subapList_dsp1[1]

		#--- LGS 0 ---
		toto=np.where(self._SH1[:,0:132].transpose() ==1)
		idx=np.array(np.where( np.in1d(toto_dsp0,toto[0]*264+toto[1]) ))[0]

		tmp=np.array(range(0,self.nsubap_dsp0*2))
		tmp0=tmp[idx*2]
		tmp1=tmp[idx*2+1]
		tmp=np.zeros(176*2)
		tmp[::2]=tmp0
		tmp[1::2]=tmp1
		remap_list.extend(np.array(tmp))

		#--- LGS 2 ---
		jump = 8 + (self.max_nb_subapPerDsp+3)*2
		toto=np.where(self._SH2[:,132:264].transpose() == 1)

		idx=np.array(np.where( np.in1d(toto_dsp1,toto[0]*264+toto[1]) ))[0]
		tmp=np.array(range(0,self.nsubap_dsp1*2))[::-1]
		tmp0=tmp[idx*2+1]
		tmp1=tmp[idx*2]
		tmp=np.zeros(176*2)
		tmp[::2]=tmp0
		tmp[1::2]=tmp1
		remap_list.extend(np.array(tmp)+jump )

		#--- LGS 1 part 1/2 ---
		toto=np.where(self._SH3[:,0:132].transpose() == 1)

		idx=np.array(np.where( np.in1d(toto_dsp0,toto[0]*264+toto[1]) ))[0]
		tmp=np.array(range(0,self.nsubap_dsp0*2))
		tmp0=tmp[idx*2]
		tmp1=tmp[idx*2+1]
		tmp=np.zeros(self.nsubap_dsp0*2-176*2)
		tmp[::2]=tmp0
		tmp[1::2]=tmp1
		remap_list.extend(np.array(tmp))

		#--- LGS 1 part 2/2 ---
		toto=np.where(np.fliplr(self._SH3[:,132:264]).transpose() == 1)

		idx=np.array(np.where( np.in1d(toto_dsp1,toto[0]*264+toto[1]) ))[0]
		tmp=np.array(range(0,self.nsubap_dsp1*2))[::-1]
		tmp0=tmp[idx*2+1]
		tmp1=tmp[idx*2]
		tmp=np.zeros(self.nsubap_dsp1*2-176*2)
		tmp[::2]=tmp0
		tmp[1::2]=tmp1
		remap_list.extend(np.array(tmp)+jump)

		init_jump=self._sc_ReplyVectorPtr+ len_hdr_RV+2*2 # first 2 dummy subapertures
		remap_list=list(np.array(remap_list)+init_jump)

		return remap_list

	def numSubapToDo(self):
		''' 
		Creates the 'dsc_NumSubapsReady', i.e. the centroid computation trigger vector,
		corresponding to the dsc_SubapPixel, previously defined.
		'''
                list_dsp0=self._subap_Numbering(self.arrayHuman[:,0:132]*1,\
                                                        list([[],[],[]]))

                list_dsp1=self._subap_Numbering(np.fliplr(self.arrayHuman[:,132:264]*1),\
                                                        list([[],[],[]]))

		numSubapToDo_dsp0=np.zeros(132,'int')
		numSubapToDo_dsp1=np.zeros(132,'int')
		for i in range(0,132-8):
			idx=np.array(np.where(np.array(list_dsp0[0]) ==[i]))
			idx1=np.array(np.where(np.array(list_dsp1[0]) ==[i]))
			if idx.any():
				numSubapToDo_dsp0[i+8]=(list_dsp0[2][idx[0][-1]])+1
			else:
				numSubapToDo_dsp0[i+8]=(numSubapToDo_dsp0[(i+8)-1])

			if idx1.any():
				numSubapToDo_dsp1[i+8]=(list_dsp1[2][idx1[0][-1]])+1
			else:
				numSubapToDo_dsp1[i+8]=(numSubapToDo_dsp1[(i+8)-1])

		numSubapToDo_dsp0=numSubapToDo_dsp0/2
		numSubapToDo_dsp1=numSubapToDo_dsp1/2

		numSubapToDo_dsp0[-1]=(self.nsubap_dsp0+3)/2
		numSubapToDo_dsp1[-1]=(self.nsubap_dsp1+3)/2

		return [numSubapToDo_dsp0,numSubapToDo_dsp1]

	def show(self):
		'''
		Show :
		- the 3 SH pupils definition
		- the Subaperture number for each DSP
		- the virtual dividing line between the 2 halves of the pnCCD
		'''
		list_dsp0=self._subap_Numbering(self.arrayHuman[:,0:132]*1,\
                                                        list([[],[],[]]))
                list_dsp1=self._subap_Numbering(np.fliplr(self.arrayHuman[:,132:264])*1,\
                                                        list([[],[],[]]))


		plt.figure()
		plt.imshow(self.arrayHuman,interpolation='nearest')

		for i in range(0,len(list_dsp0[0])):
			plt.text(list_dsp0[0][i],list_dsp0[1][i], str(list_dsp0[2][i]),\
					 color='w',verticalalignment='top',\
					 size=8)
		
		for i in range(0,len(list_dsp1[0])):
			plt.text(264-8-list_dsp1[0][i],list_dsp1[1][i], str(list_dsp1[2][i]),\
					 color='w',verticalalignment='top',\
					 size=8)
		
		plt.axvline(132, color='green', linewidth=4)
		plt.title('BCU subaperture definition, DSP0 &DSP1')

		plt.show()



	def _createSubapDef(self, centers, radii):
		'''
		Enable to create the 3 SH puils definition.

		center=[c1,c2,c3]
		radius=[r1,r2,r3]
		- c1, c2, c3 are [x,y] central coordinates of SH pupils;
		- r1, r2, r3, radius in # of subap per diameter

		  1: left bottom, 2: right bottom, 3: top

		 also need a check of subap number < max per dsp;
		 and no conflict between 2 DSPs

		 typically: r=14 here (offset of one for code reason)
		'''
		c1=centers[0]
		c2=centers[1]
		c3=centers[2]
		r1=radii[0]-1
		r2=radii[1]-1
		r3=radii[2]-1
		array=np.zeros((264,264))
		array_sh1=np.zeros((264,264))
		array_sh2=np.zeros((264,264))
		array_sh3=np.zeros((264,264))
		array2_sh1=np.zeros((264,264))
		array2_sh2=np.zeros((264,264))
		array2_sh3=np.zeros((264,264))
		idx=np.indices((264,264))
		sh1=np.sqrt((idx[0]-c1[0])**2 + (idx[1]-c1[1])**2)
		sh2=np.sqrt((idx[0]-c2[0])**2 + (idx[1]-c2[1])**2)
		sh3=np.sqrt((idx[0]-c3[0])**2 + (idx[1]-c3[1])**2)

	        array_sh1[c1[0]-4-r1*4:c1[0]+4+r1*4,c1[1]-4-r1*4:c1[1]+4+r1*4]=1
	        array_sh2[c2[0]-4-r2*4:c2[0]+4+r2*4,c2[1]-4-r2*4:c2[1]+4+r2*4]=1
	        array_sh3[c3[0]-4-r3*4:c3[0]+4+r3*4,c3[1]-4-r3*4:c3[1]+4+r3*4]=1

		# SH top left
		tmp=array_sh1*1
                idx=np.where(tmp==1)
                x_ref= idx[0][0]
 		y_ref=idx[1][0]

		while np.sum(tmp)>0:
			x_ref= idx[0][0]
			y_ref=idx[1][0]
		
			sub_radii=sh1[x_ref:x_ref+8,y_ref:y_ref+8]
			if np.mean(sub_radii)  > r1*4+4 :
				array_sh1[x_ref:x_ref+8,y_ref:y_ref+8]=0
		
			tmp[x_ref:x_ref+8,y_ref:y_ref+8]=0
			idx=np.where(tmp ==1)


		# SH top right
		tmp=array_sh2*1
		idx=np.where(tmp ==1)
                x_ref= idx[0][0]
 		y_ref=idx[1][0]
		while np.sum(tmp)>0:
			x_ref= idx[0][0]
			y_ref=idx[1][0]
			sub_radii=sh2[x_ref:x_ref+8,y_ref:y_ref+8]
			if np.mean(sub_radii)  > r2*4+4 :
				array_sh2[x_ref:x_ref+8,y_ref:y_ref+8]=0
			
			tmp[x_ref:x_ref+8,y_ref:y_ref+8]=0
			idx=np.where(tmp ==1)

		# SH bottom middle
		tmp=array_sh3*1
		idx=np.where(tmp ==1)
                x_ref= idx[0][0]
 		y_ref=idx[1][0]
		while np.sum(tmp)>0:
			x_ref= idx[0][0]
			y_ref=idx[1][0]
			
			sub_radii=sh3[x_ref:x_ref+8,y_ref:y_ref+8]
			if np.mean(sub_radii)  > r3*4+4 :
				array_sh3[x_ref:x_ref+8,y_ref:y_ref+8]=0
			
			tmp[x_ref:x_ref+8,y_ref:y_ref+8]=0
			idx=np.where(tmp ==1)


		# central subaperture is vignetted
		array_sh1[c1[0]-4:c1[0]+4,c1[1]-4:c1[1]+4]=0
		array_sh2[c2[0]-4:c2[0]+4,c2[1]-4:c2[1]+4]=0
		array_sh3[c3[0]-4:c3[0]+4,c3[1]-4:c3[1]+4]=0
		
		# Adding together the three SH pupils
		array=array_sh1 + array_sh2 + array_sh3

		# Updating internal variable of SH1, SH2, SH3
		self._SH1 = array_sh1
		self._SH2 = array_sh2
		self._SH3 = array_sh3

		if np.where(array > 1)[1].size !=0:
			warnings.warn("Wrong pupil definitions: Shack-Hartmann are overlapping")
			print "Wrong pupil definitions: Shack-Hartmann are overlapping"

		return array


	def _array2subap_pix_ptr(self,tmp):
		'''
		converts an 264x264 array with pupil definition (px=1, =0 elsewhere) to the dsc_SubapPixel
		format as requested for the BCU
		'''
		array_sh_coord_dsp0=tmp[:,0:132]*self._arrayIndices()
		array_sh_coord_dsp1=np.fliplr(tmp[:,132:264])*self._arrayIndices()

		dsp0=self._recursive_defSubapBCU(array_sh_coord_dsp0,list([]),\
							 1,self.nsubap_dsp0,0,self.subapSize)
		dsp1=self._recursive_defSubapBCU(array_sh_coord_dsp1,list([]),\
							 1,self.nsubap_dsp1,1,self.subapSize)
		
		dsp0=list(np.array(dsp0,int))
		dsp1=list(np.array(dsp1,int))
		# copy 1st, 2nd and last subaperture... see manual

		header0=dsp0[:128]
		header1=dsp1[:2*64]

		dsp0.reverse()
		header0.reverse()
		dsp0.extend(header0)
		dsp0.reverse()

		dsp1.reverse()
		header1.reverse()
		dsp1.extend(header1)
		dsp1.reverse()

		if np.mod(self.nsubap_dsp0,2) :
			footer0=dsp0[-64:]
			tmp=np.zeros(128,'int')
			tmp[::2]=footer0
			tmp[1::2]=footer0
			footer0=tmp
			dsp0[-64:]=footer0[:64]
			dsp0.extend(footer0[64:128])

		if np.mod(self.nsubap_dsp1,2) :
			footer1=dsp1[-64:]
			tmp[::2]=footer1
			tmp[1::2]=footer1
			footer1=tmp
			dsp1[-64:]=footer1[:64]
			dsp1.extend(footer1[64:128])

		vectorsBcu=[dsp0,dsp1]
		return vectorsBcu

	def _arrayIndices(self):
		'''
		we assume 264x264 array pnCCD format.
		'''
		array_dsp=np.zeros((264,132))

		array_dsp[0:66,]=np.array(range(0,264*132,4)).reshape(132,66).transpose()
		array_dsp[66:132,]=np.array(range(1,264*132,4)).reshape(132,66).transpose()
		array_dsp[132:198,]=np.array(range(2,264*132,4)).reshape(132,66).transpose()
		array_dsp[198:264,]=np.array(range(3,264*132,4)).reshape(132,66).transpose()

		return array_dsp
	

	def _recursive_defSubapBCU(self,array_coord_shdef,subap_pix_vector,counter,\
					   lastSubnb,dsp,size):
	        '''
		used by _array2subap_pix_ptr. Enable to recursively find and define pixel of a subaperture.
		Do it '2by2' subaperture (as requested for the BCU for optimization.
		'''
	        #counter+=2
		array_coord_shdef=array_coord_shdef*1 

		assert 34848 == np.size(array_coord_shdef)
	
		array_coord_shdef[0,0]=-1
		x_ref=np.where(array_coord_shdef.transpose()>0)[0]
	
		if (not(len(x_ref)>0) or counter>275):
			pass
		else:
			x_ref=np.min(x_ref)
			y_ref=np.min(np.where(array_coord_shdef.transpose()[x_ref,]>0))
	
			tmp=(array_coord_shdef.transpose()[x_ref:x_ref+size,\
								  y_ref:y_ref+size])[::(-2*dsp+1),:]
			#print tmp
                        tmp0=(tmp.flatten())

			array_coord_shdef.transpose()[x_ref:x_ref+size,y_ref:y_ref+size]=0	
			counter+=1		
			# *** defining second subaperture
			if np.mod(lastSubnb,2) and counter==(lastSubnb+1):
				tmp2=tmp0
			else:
				x_ref=np.where(array_coord_shdef.transpose()>0)[0]
				x_ref=np.min(x_ref)
				y_ref=np.min(np.where(array_coord_shdef.transpose()[x_ref,]>0))
				tmp=(array_coord_shdef.transpose()\
				    [x_ref:x_ref+size,y_ref:y_ref+size])[::(-2*dsp+1),:]
                        	tmp1=(tmp.flatten())
				tmp2=np.zeros(128)
				tmp2[::2]=tmp0
				tmp2[1::2]=tmp1

			array_coord_shdef.transpose()[x_ref:x_ref+size,y_ref:y_ref+size]=0			
			counter+=1
			subap_pix_vector.extend(list(tmp2))	
			self._recursive_defSubapBCU(array_coord_shdef,subap_pix_vector,counter,lastSubnb,dsp,size)
	
		return subap_pix_vector



	def _subap_pix_ptr2array(self,m_subap_pix_ptr_dsp0,m_subap_pix_ptr_dsp1, values=[[1],[1]]):
		'''
		Transform an subap_pixel_ptr as defined for the BCU into a 
		normal array (for 'human representation')
		'''
		
		m_subap_pix_ptr_dsp0=np.array(m_subap_pix_ptr_dsp0)
		m_subap_pix_ptr_dsp1=np.array(m_subap_pix_ptr_dsp1)
		
		pixelvec0=np.zeros(264*132)
		pixelvec0[m_subap_pix_ptr_dsp0]=values[0]
		
		pixelvec1=np.zeros(264*132)
		pixelvec1[m_subap_pix_ptr_dsp1]=values[1]
		
		array=np.zeros((264,264))
		for i in range(0,132):
		    array[0:66,i]=pixelvec0[(0+i*264):(261+i*264):4]
		    array[66:132,i]=pixelvec0[(1+i*264):(262+i*264):4]
		    array[132:198,i]=pixelvec0[(2+i*264):(263+i*264):4]
		    array[198:264,i]=pixelvec0[(3+i*264):(264+i*264):4]
		
		c=263
		for i in range(0,132):
		    array[0:66,c-i]=pixelvec1[(0+i*264):(261+i*264):4]
		    array[66:132,c-i]=pixelvec1[(1+i*264):(262+i*264):4]
		    array[132:198,c-i]=pixelvec1[(2+i*264):(263+i*264):4]
		    array[198:264,c-i]=pixelvec1[(3+i*264):(264+i*264):4]
		
		return array
		

	def _subap_Numbering(self,array,subap_number_list,counter=0,size=8):
		'''
		used by numSubapToDo, and remapVector and remapVectorBCU
		'''
		array[0,0]=-1
		x_ref=np.where(array.transpose()>0)[0]
	
		if (not(len(x_ref)>0) or counter>275):
			pass
		else:
			x_ref=np.min(x_ref)
			y_ref=np.min(np.where(array.transpose()[x_ref,]>0))
	
			subap_number_list[0].extend(list([x_ref]))
			subap_number_list[1].extend(list([y_ref]))
			subap_number_list[2].extend(list([counter]))
		
			array.transpose()[x_ref:x_ref+size,y_ref:y_ref+size]=0
			counter+=1
			self._subap_Numbering(array,subap_number_list,counter)
	
		return subap_number_list

	

if __name__ == "__main__":

    subap=SubapDef()
    print "OK"		
		
		
		
	
